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We investigate the thermal entanglement of interacting two qubits. We maximize it by tuning 
a local Hamiltonian under a given interaction Hamiltonian. We prove that the optimizing local 
Hamiltonian takes a simple form which dose not depend on the temperature and that the corre- 
sponding optimized thermal entanglement decays as 1/{T log T) at high temperatures. We also find 
that at low temperatures the thermal entanglement is maximum without any local Hamiltonians 
and that the second derivative of the maximized thermal entanglement changes discontinuously at 
the boundary between the high- and low-temperature phases. 

O ! I. INTRODUCTION 
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Quantum entanglement plays an essential role in quantum information processing [l[ . Various kinds of investiga- 
tion have been carried out to understand properties of entanglement for the last two decades 0, Q. The thermal 
entanglement which is entanglement of thermal equilibrium states, is one of the important concepts because it 
' shows us the effect of thermal fluctuations on entanglement. Thermal disturbances generally cause disentanglement 
and have serious effects on quantum information processing. Therefore, many schemes have been proposed to protect 
(— I ' entanglement from thermal disturbances |5l-[l4|. As one of these schemes, a lot of attention has been paid to methods 
O ^. based on manipulation of local Hamiltonians [sl [sl. [Tol - [l^ : for example, in quantum spin systems, bipartite thermal 
^ I , ' entanglement can be enhanced by modulating external magnetic fields. In the present paper, we focus on a simple 
P5 \ question as to how much entanglement can be generated by optimizing the local Hamiltonian. We give a theoretical 

■ limit of entanglement enhancement by manipulation of the local Hamiltonians. 
^ Relationships between the thermal entanglement and local parameters have been investigated especially in bipartite 

quantum spin systems From these researches, behavior of the thermal entanglement under external 

magnetic fields may be understood in the cases of almost all interactions. However, little has been reported on the 
maximization problem of the thermal entanglement; in the case of the bipartite XY spin model, this problem has been 
^ ' solved only numerically 12 1. Until now, there are no analytical approaches to optimizing the thermal entanglement 



00 



o 



of arbitrarily interacting two qubits. 

In the present paper, we will answer the following question: given a system of two qubits which interact via an 
[ arbitrary interaction Hamiltonian, how can wc maximize the thermal entanglement between these two qubits by 
changing only the local Hamiltonian? A naive approach to this problem may be to solve the optimization problem 
numerically. However, this problem has six local parameters in total and the functional forms of entanglement measures 
such as the concurrence (l9| and the negativity (20| are very complicated. Thus, for an arbitrary interaction, it is 
^-H ' difficult to solve this optimization problem numerically. Therefore, we employ perturbation techniques and utilize 
r* , symmetric properties in order to determine the optimizing local Hamiltonian analytically. In this way, for all kinds 
. !^ of interaction, we give general properties of the optimized entanglement. 
] Our main results are the following: 
H ■ 

I 1. We find that at low temperatures the thermal entanglement is maximum without any local Hamiltonians, 

whereas at high temperatures it is maximized by non-zero local fields. We refer to the former temperature range 
as the low-temperature phase and the latter temperature range as the high-temperature phase. The secondary 
differentiation of the maximized entanglement is discontinuous at the phase boundary. 

2. In the high-temperature phase, the functional form of the optimizing local Hamiltonian is independent of the 
temperature; only the coefficients depend on the temperature. 

3. The optimized entanglement, enhanced by a local Hamiltonian in the high-temperature phase, decreases with 
increasing temperature as l/(TlogT). 

4. If the interaction Hamiltonian has no degeneracy of its eigenvalues, the entanglement is maximized without local 
Hamiltonians over a finite range of the low-temperature phase. 

5. If the interaction Hamiltonian has degeneracy, the low-temperature phase shrinks to the zero-temperature point. 
The optimizing local Hamiltonian becomes infinitesimal and the optimized entanglement becomes full in the 
low-temperature limit. 
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Our paper is organized as follows. In Section II, we state the main problem after symmetry consideration. In 
Section III, we give the main theorems on the entanglement optimization. In Section IV, we show numerical results 
of the optimizing local parameters, the boundary temperatures and the singularity at the phase boundary. We also 
argue that the two phases appear because of competition between the purifying effect and the decoupling effect both 
of the local Hamiltonian. Finally, in Section V, a discussion concludes the paper. 



II. ENTANGLEMENT OPTIMIZATION PROBLEM 



First, we set the fundamental framework of the present problem. We consider a 2 (g) 2 system of cti and 02- The 
most general form of the Hamiltonian of this system is given as follows: 

i?tot = ^int + ^^LO 7 

i^LO^ Yl (hWi^I + ht^I^a'^), (1) 

i—x,y,z 

where {(y\}i=x.y.z and {<^2\i=x,y,z are the Pauli matrices, i?int is an interaction Hamiltonian, and i?LO is a local 
Hamiltonian. We assume that {Jij}ij=x,y,z are fixed and independent of the temperature, whereas we can change 
the parameters {hf, h\, /if, /if, /I2, /ill arbitrarily. 

We parametrize the local fields in the polar coordinates as 

{h\}i=x.y.z ~ {hi sin^i cos 01, hi sinOi sin 01, hi cos^i}, 
{h2}i=x.y.z = {h2 sin 6*2 COS02, /i2 sin 02 sin02, /i2 cos 6*2}. 

(2) 



Hereafter, we use the parametrization 



, hi + h2 . hi - h2 

h^ o ' C = , . , , (3) 

2 hi+ h2 



where — 1 < C !i 1 and /i > 0; in other words, 

hi = {l + C)h, h2 = {l~0h. (4) 

Then, the four eigenvalues of ffLO are 

{-2/i,-2C/i,2C/i,2/i}, (5) 

where we define the corresponding eigenstates as {| ), | — h), |H — ), 1++)}- 

The density matrix in thermal equilibrium is 

P=^, (6) 

where Z = tr(e~^^'°') is the partition function and /3 = 1/kT with k the Boltzmann constant. In order to quantify 
entanglement, we adopt the negativity j20j as an entanglement measure. The negativity is defined as the trace norm 
of a partially transposed density matrix: 

N{p) = \\p'^^\\i-1 

= max(-2A_,0), (7) 

where || ||i is the trace norm, Ti denotes the transpose with respect to only ci, and A_ is the minimum, possibly 
negative eigenvalue of p^^ . The second equation of ([7]) comes from the fact that p^^ can have only one negative 
eigenvalue, if any (2]| . Thus, the present entanglement optimization problem is equivalent to finding the values of 
{hf, h\, hi, /if, /12, /12} which maximize N(p) for an arbitrary fixed interaction Hint- 

Before presenting our main results on the entanglement optimization, we prove the following Lemma 1 to simplify 
the present entanglement optimization problem. 
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Lemma 1. By local unitary transformations of i?int, we can eliminate the interaction parameters {Jij}i^j and 
reduce it to the form 

Hint= '^^'^i®^- (8) 

i—x,y,z 

We can also choose the parameters {Jx, Jy, Jz} such that {J^^Jy} > Jz > or > Jz > {Jx,Jy}- 

In spin-1/2 systems, this means that we can transform any interactions including the Dzyaloshinskii-Moriya 
(DM) jl4l I22I |23| interaction into a ferromagnetic or an anti-ferromagnetic Heisenberg exchange interaction. 

Proof. We can prove this Lemma by applying a singular value decomposition [23 | to the matrix {J)ij = Jij. In 
this case, the singular value decomposition UJW is performed by 3 x 3 real orthogonal transformations U and W of 
the three-dimensional spin spaces of the spins 1 and 2, respectively. A real orthogonal transformation is composed of 
rotation and inversion operations, but inversion operations cannot be performed by unitary transformations. There- 
fore, we remove the inversion operations from the real orthogonal transformation of the singular value decomposition 
and restrict ourselves only to the rotation operations, which means detU = detW =1. In other words, we rotate 
(Ti — {erf , erf , erf} with U and a2 = {"'f j "'2 j "'2} with W. Then we can transform {Jij}ij=x,y,z into the antiferromag- 
netic cases {J^, Jy} > J2 > or the ferromagnetic cases 0> Jz > {Jx, Jy}, with the other elements {Jij}i^j put to 
zero. Here, we choose the z-axis so that | J^l is the least of {\Ji\}i=x,y,z- Thus, Lemma 1 is proved. 

Let us show an example in the case of the XXZ model with the z-component of the DM interaction. The 
Hamiltonian of such a system is given by 

^^int = Jcrf ® (T2 + <^ cr^ + JzCrt 

+ Dz{Gl®Gl-a\®al), (9) 

where J and Jz are the real coupling cocfhcicnts and Dz is the z-component of the DM interaction. In the case of 
J = 1, Jz = —2 and Dz = 1, we can transform {Jxi Jy^ Jz, Dz} into { — v^, — v^, —2, 0} by rotating the spin 1 by 135 
degrees around the z-axis, namely into 

Hi-^t = - V2al ® erf - V2erf ® cr| - 2crf ® (72 ■ (10) 

This is an antiferromagnetic Heisenberg interaction. To attain this result, first, the singular value decomposition 
transforms { J^, Jy, Jz,Dz} into {V^, V^, 2, 0} by rotating the spin 1 by —45 degrees around the z-axis and inverting 
the z-axis of the spin. Next, we remove the inversion of the z-axis because it cannot be performed by unitary 
operations, and thereby transform {J^, Jy, Jz, Dz} into {V^, \/2, —2, 0}. By changing the rotation angle from —45 to 
135, we can invert the signs of J^ and Jy and arrive at { J^, Jy, Jz, Dz} = {— V2, — V2, —2, 0}. 

In the following, based on Lemma 1, we always use the diagonalized form ([5]) of the interaction parameters with 
{Jx,Jy} > Jz > or > Jz > {Jx,Jy}- We now have all the necessary ingredients to state the main theorems. 



III. MAIN ANALYTICAL RESULTS 



In the present section, we analytically discuss the optimization problem. The main conclusion of the present section 
is that the negativity is maximized by the parameters {hf, h\, /if, /if, = {0, 0, ft-op, 0, 0, —hop}- The optimizing 
parameter ft,op must be very large at high temperatures, whereas it may be at low temperatures. 



A. Optimization in the high-temperature limit 

Let us first discuss the optimization problem in the high-temperature limit. 

Theorem 1. In the high-temperature limit /3 — > 0, the local parameters which maximize the entanglement N[p) are 
given in the form of {/if, h\, h\, /ij , h^, /12} = {0> 0, /lop, 0, 0, — /lop}- The optimizing value /lop is given by the solution 
of the following equation: 

where 

Kp = /3/iop, (12) 
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and the optimized entanglement A'op asymptotically behaves as 

\Jx 

op op 



where we used Eq. upon moving from the first line to the second line. 
The leading order of the solution of Eq. ((TT|) is given by 



We can thereby obtain the following simpler asymptotes: 

/lop^i^l^ as ^^0, (15) 

^op(p)^/3^t^ as P^O. (16) 

That is, the optimizing value hop depends only on the temperature and the optimized negativity decays in the form 
l/(riogr) in the limit /3 — >■ 0. In Appendix A, we compare the asymptotes of Eqs. ([TT|) and with those of 
Eqs. (113 and ^ 

Proof. We prove Theorem 1 in the following steps. First, we prove in Lemma 2 that the optimizing local parameter 
hop is greater than or equal to (log 1/ /3)/{2f3) in the high-temperature limit and the optimized thermal state is nearly 
a pure state. The entanglement of the state comes from perturbations to the pure state. Then, we calculate the neg- 
ativity approximately by perturbation method in Lemma 4. Using this expression, we finally solve the maximization 
problem for each local parameter. 

First, we determine a lower bound of the optimizing value hop and prove that the optimized thermal state is a 
nearly pure state. For this purpose, we prove the following Lemma 2. 

Lemma 2. A necessary condition for the existence of the entanglement in the high-temperature limit under a fixed 
interaction Hamiltonian Hi^t is given by 

j3h > as /3 ^ 0. (17) 

This Lemma 2 shows that (log l//3)/(2/3) is a lower bound of the optimizing value of hop- 
Proof. We firstly prove that we need a non-zero value of /3h for the existence of the entanglement in the high- 
temperature limit /? — > 0. In other words, we need h at least of order 1/0. In order to show this, we consider a general 
necessary condition for the existence of the entanglement given by |25| 



Ai > A3 + 2VA^ > 3A4, (18) 

where {A^}^^]^ are the eigenvalues of the density matrix p in the non-ascending order (Ai > A2 > A3 > A4). Let us 
define the eigenvalues of Htot = -f^LO + ^^int as in the non-descending order {Ei < E2 < E3 < E4). Equa- 

tion ^ gives the eigenvalues of p as {e^^^^' /Z}'^^^, and therefore the inequality p^ . or e^^^^ > Se"^-^*", gives 

/3iEi-Ei)>\og3. (19) 

Here, Hi^t is a constant matrix and hence /3-ffint — !• as /3 0. If we let i?LO be of the same order as Hint, the 
left-hand side of p9|) would vanish in the limit /3 and (fT8|) would not be satisfied. Therefore, we have to make 
Hlo much greater than i^int, and then the eigenvalues of Htot should converge to those of Hi^o, {— 2/i, —2(^h, 2(^h, 2h} 
in the limit /? — > 0. With Ei — > —2h and E4 — > 2h, the inequality (fT9|) reduces to the following inequality: 



Ph>'-^. (20) 



This inequality means that we need a non-zero value of I3h in the high-temperature limit /? — )■ 0. In other words, we 
need to make h grow as 1//3 at least, in order for the entanglement to exist in the limit /3 — > 0. 
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Next, we derive an approximation of the density matrix, and then obtain Eq. (|17p by utihzing the Peres- Horodecki 
criterion (26l. [27|. which is a necessary and sufficient condition for the existence of the entanglement. In the present 
optimization problem, we fix Hint to a constant matrix, and therefore we have /3i?int — > in the high-temperature 
limit. We thereby work in the first-order approximation with respect to ^Hi^t' 

Zp = e--f^Lo-/3Hint 
Jo 

= e-''-o ^ pJ2UA^^\HintW)\^^)W\, (21) 



where Z = tr(e 



), and we let H{^q = (3Hi^o with h' = /3h as well as 



-E' 



for E'^^El, 
iovE'^=El. 



(22) 



Here, {E'^y^^^ are the eigenvalues of H!^q = /Si^LO, {—2h',~2C,h',2C,h',2h'}, and are the corresponding 

eigenstates, {| ), | — h), !++)}• 

We then utilize the necessary and sufficient condition for the existence of the entanglement, det p^^ < 0. This 
has been proved [2l| to be equivalent to the Peres- Horodecki criterion [2^, (2^ . In the following discussion, among 
the various terms of the expansion of det p'^^ , we compare the values of the products including off-diagonal elements 
(POD) with that of the product of the diagonal elements (PD), which has a positive value. Then a necessary condition 
for det < is that POD is greater than or of the same order as the PD. 

To analyze the order of the PD and the PODs, we express p^^ in the basis {\p)}fi=i and focus on the main terms 
for C 7^ and C ^ ±1: 



Zp 



"21 — 



023 



"31-ft7- 
n 

o-^x—jT- 

P-2C/1' 



032 — 



042 
O34 



fie 



2Qh' 



-2C/i' 



O43 



Pe 



-2Ch' 



h' 
-2h' 



(23) 



where {o,j} are determined from Eqs. ([?!]) and and are constants of order 1. Note that on the diagonal of 
Eq. ()23p . the second term in Eq. ()2ip is neglected in comparison to the first term. Then we compare the orders of 
the PODs with that of the PD. The PD is given by g'^^' ^2C,h' ^-'2C,h' ^-2h' ^ whereas each POD includes at least two 
off-diagonal elements. The maximum of the absolute value of the PODs is of order e*'' /3^//i'^, which comes from the 



product — e^'* x 041 x 014 '^'^, x e . Therefore, it is necessary for dot p'^^ < that e""* f3'^/h''^ is greater or of 

order 1, which is the order of PD. By taking the logarithm of e""* we can obtain the following inequality as a 

necessary condition: 



> 



> 



> 



logl//3 logfe^ 
2 2 

logl//3 
2 

logl//3 
2 

logl//3 



1 Zlogl//? logV 
- log 

2 2 2 



ilogfi^^ + ilogfi^^ 
2 2 2 4 



(24) 



where we utilized (|20p in deriving the third inequality and used the fact /? — >■ in deriving the last inequality. Thus, 
Lemma 2 is proved for C 7^ and C 7^ ±1- For C = OorC = ±l, some of the eigenvalues of H^q are degenerate, which 
means that E'^ can be equal to E^ in Eq. ([22)) . and Zp'^^ is not of the same form as that of Eq. ([23]). However, the 
inequality ([M|) still holds as is proved in Appendix B. 

We now consider the negativity ([7]) in the range given by (IT71) . We first show in the following Lemma 3 that the 
optimized negativity in the cases of ^ = ±1 is not large enough. 
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Lemma 3. In the cases of C 



±1, the optimized negativity satisfies the fohowing: 

iVop(p,c = ±i) 0- 



/3 



^ 0. 



(25) 



This lemma shows that the optimized negativity in the cases of C = ±1 is of a higher order of /3. Indeed, we numerically 
confirmed in the cases of C = ±1 that the entanglement exists, but its amplitude is of order /3^. 

Proof. Let us prove Eq. ([25]) in the case of C = 1- The proof for C = ^ 1 is almost the same. We start from the 
main term of Z p^^ for C = 1 in the representation in the basis 



Zp 



( e^'*' ai2/3e2''' 



ai3 



h' 



an — 

«24 — 



032 \ 



a43/3e 



-2h' 



a42-^ 
034/36-2" 
^-2h' 



(26) 



where we used the fact that at C = 1 the eigenvalues of H{^q in Eq. ([22]) are degenerate as {E[, E'2, E';^, E'^} = 
{2h',2h',~2h',-2h'}. 

In order to optimize the negativity, we necessarily consider the region h' = /3h > (log l//3)/2 as is given in Lemma 2. 



Therefore, we can use the fact e^'* > I3~^ in (|26|) . Of the elements of the matrix ([26)), the (1, 1) and (2, 2) elements 
are of order /3~^ or greater, whereas the (3,3), (4,4), (3,4) and (4,3) elements are of order f3^ or less. The other 
elements are approximately of order 1. We therefore break up the matrix (|26p in the form 



1 
Z 





/^2h' 


o\ 






1 





e^'^' 






z 




I 



0^ 






( 







pe^'^' 


032 \, 







pe^'^ 
041 ^, 


Pe-^" 
042 h' 




h' 










Pe^" 





/ 



(27) 



where Z ~ 2e^^ + 2e ^'^ ~ 2//3, and therefore the first term is the dominant term of order 1, whereas the second 
term is of order The eigenvalues of the dominant term are {e^^ /Z, e^^ /Z, 0, 0} and the corresponding eigenstates 

are {| ), | — h), |H — ), |++)}. A negative eigenvalue can appear when the degeneracy of the two zero eigenvalues of 

the states |H — ) and |H — h) is resolved by perturbation. Then, the level repulsion between them makes one of them 
positive and the other negative. However, the first-order perturbation of the second term of Eq. (P7)) dose not resolve 
the degeneracy of the zero eigenvalues. Therefore, the negative eigenvalue must be produced in a higher order of /? in 
the case of C = 1. Thus, Lemma 3 is proved. We focus on the cases ^ 7^ ±1 hereafter. 

Using the lower bound (|17|) of the optimizing parameter ft,op, we next prove that the optimized thermal state 
is a nearly pure state in the cases of ^ 7^ ±1. For this purpose, we consider the eigenstates of the perturbed 

density matrix. We define the perturbed eigenstates of fiHtot = H^o P^^^t as {| '),| — +'),|H — 

corresponding to the eigenstates {| ),| — +),|H — ),!++)} of ^lc respectively, and their eigenvalues as 

{2h' — (36ei,2C,h' — pSe2,—2(h' — f^Se^, —2h' — /3(Se4}, where {(5ei}i=i are the perturbative changes due to Hint, which 
are of order 1. Then the density matrix is given by the summation over these four states. In the high-temperature 
limit /? — 7> 0, the mixing ratio {A__',A /,A++'} of the states {| '), | — h'), \ +—'), |++')} are 



{A — ', A_+', A+_/, A++'} 
1 

Z' 



^_^p2h'-pSei ^ g2C/i'-/3<5c2 ^ ^-2Ch'-p5e3 ^ g-2/i'-/J5£4 | ^ 



(28) 



where 



2 ^ ^2h'-pSei ^2Ch'-pSe2 _^ ^-2Ch'-pS€3 ^ ^-2h'-pSe 



(29) 
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In the region h' > (log l//3)/2, which is the lower bound of h' , we have 



A- 



A- 



A— 



^ g-2(l-C)'»'+/3<5ei 



^ g-2(l+C)h'+/35ei- 



-l3Se2 



</3 



l-Cg-/3<5«i+/3<5«2 



-2/i'+;35ei-^(5c4 



(30) 



Since the right-hand sides of the inequalities vanish in the limit /3 — >■ 0, we deduce that the optimized thermal state 

is a nearly pure state of | ') in the high-temperature limit /3 — when ^ ^ ±1. 

Next, we perturbatively calculate the negativity in the cases of C 7^ ±1. Since the optimized state is a nearly pure 
state of I — 



we regard the other contributions { | 

Po = I 
6p = 



-')(- 



E 

{i=-+', + -'.++'} 



-+')} as perturbation: 
Kpi, 



(31) 



where {| '),| — h^ilH — the eigenstates of /3i?tot ^ H^q + l3Hint as has been stated. In order to 

calculate the negativity approximately, we derive the expression for the perturbation of the negativity caused by an 
infinitesimal variation of the density matrix. 

Lemma 4- When the negativity has a non-zero value, the first-order perturbation of the negativity is given by 



N{pa + Sp) ~ N{po) 



-\Sp 



(32) 



where we refer to the eigcnstate corresponding to the negative eigenvalue of p^^ as 

Proof. The non-zero negativity is given by the negative eigenvalue A_ of the partial transpose of the density matrix, 
Pq^ , as is defined in ([7|). Because of the linearity of the partial transpose, if po changes into po + Sp, Pq^ also changes 
into Pq^ + dp^^ . Moreover, the eigenstate of p^^ corresponding to A_ is not degenerate because A_ is the only possible 
negative eigenvalue [2l|. Then, from the general perturbation theory for A_, we have Eq. p2p in the first order. 

From Eqs. (PT|) and ([5^ . we can calculate the negativity in the present case of ^ 7^ ±1 as 



n(J2Kp.) 

'))- 



-N{\ 



E 



2A,; 



6_|pfM'^-)+0(A,2))- 



(33) 



The state | ') and its negativity A^(| ')) are calculated in the first order of the perturbation H{^q — H{^q + /3-ffint- 

The zeroth-order eigenstates and eigenvalues are {| ), | — +), |H — ), |++)} and {—2h' — 2(h', 2(h', 2h'}, respectively. 

The first-order eigenstate for the state | ) is then given by 



I—') = I—) + I3n,\-+) + I3n2\+-) + I3ns\++) + 0{l3^), 



(34) 



where 



111 



«2 



"3 



-2(1-C)/J' ' 

(+-|gint|— ) 

-2{C+l)h' ' 

(++|gint|— ) 

-Ah' 



Note that the normalization factor of the state 
is therefore given in the basis of {| ), | — h), 



isl + 0(/3^). The matrix representation of = (|- 
1++)} as follows by ignoring the terms of 0(/3^): 



-')(— 1) 



(35) 



(I— 'X— 1) 



Ti 















1" 




"2 

















ni 





"3 




(= 










n* 


n* 















0/ 




u 








0/ 



(36) 
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The zeroth-order eigenvalues of p^^ are {1,0,0,0}. The negative eigenvalue emerges when the degeneracy of the 
first and second zero eigenvalues resolve in the first order of (3. The third zero eigenvalue remains to be zero. The 
eigenvalues are then given by {1, /3|n3|, — /Sjrial, 0} in the first order and hence the negative eigenvalue — /3|n3| gives 
the negativity 

Ni\—'))=2/3\n,\ 



=/3 '^++'f;r'— (37) 



The corresponding eigenstate is given by 



Similarly, we have 



as well as 



A-+' 



l-^-) = -^{1-+) ^1 + -)) + Om—)- (38) 



P% = (|-+')(-+'|)'^^ = |-+)(-+| + 0(/3), 

p?-' = (l+-')(+-'l)'^^ = !+-)(+-! + 0(/3), 

P+V = {\++'){++'\)'^' = |++)(++| +0(/3), (39) 

„2C?i'-/3<5e2 



z 

= e-2(i-C)'»'+0(^2-C-|CI)^ 

2Ch'-l3S€3 

A+_, = ~ e-^i^+Oh' ^ o(^2+c-ICI)^ 

-2/i'-^(5e4 

A++' = — ^ ^e-^'''+0(/33-ICI), (40) 

where we used Eq. (p9)) for Z. Note that the first term of each of A ' and A++' is of order /J^'^^ and 

or less, respectively, in the range of PT| . h! > (logl//3)/2. By substituting Eqs. ([57)) -(|in )) in Eq. we have 



1=1 



2h' 

e-2(i-C)/>' _ e-2(i+C)'.' + 0(/32(i^|CI)) (41) 



for C 7^ ±1. 

Because the matrix element |(H — h|i?int| )| is independent of h' and C, we can solve the maximization problem of 

Eq. (j4ip as follows. First, to maximize the negative terms in Eq. (|4ip . we must put C — 0. Then, by differentiating 
Eq. ()4ip with /i', we have the optimizing parameter h'^^ as a solution of 



Q/,/2 

e^K^ f42) 

- P\{++\H,^,\—)\- ^ ' 



The optimized negativity is then given by 



where we used Eq. (|42p upon moving from the first line to the second line. This is the result for C 7^ ±1. From 
Lemma 3, we see that the optimized negativity in the case of C — is larger in the limit /3 than the one ([^S)) 
in the cases of C = ±1. 
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The other optimizing parameters to be fixed are {6*1, 0i, 6*2, 02}- Let us see how these parameters affect the value 

of P5)) . These parameters affect the matrix element |(++|iJint| )| and hence the value of P5)) directly as well 

as indirectly through h'^^ given by Eq. (|42|) . We can write down the solution of Eq. (j42|) in terms of Lambert's W 
function 1281 1 . which is defined solution of 



because we can cast Eq. (|42l) into the form 
The appropriate solution of Eq. (|^ is given by 



/3|(- 



hL 
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(44) 



(45) 



(46) 



where W-i{x) is the branch of W{x) satisfying W-i{x) < —1 in the domain — 1/e < x < [2^. The function 
—W-i{—x) is a monotonically decreasing function of x in the domain < a: < 1/e. Therefore, maximizing the element 

|(H — hl-ffintl )| with respect to the parameters {^i, 0i, ^2j ^2} brings h'^^ to its minimum within the range h'^^ > 1. 

Since the factor 



1 



2h' 



1 



1 



2h' 



(47) 



in Eq. (j43p is a decreasing function of h'^^ for h' > 1, minimizing /i^p within the range h'^^ > 1 brings the factor (j47p 

to its maximum. To summarize, the element |(h — h|-ffint| )| increases the value of (|43l) not only directly but also 

through /Iqp indirectly. 

The next task is then to find the parameters {^1,01,^2,02} that maximize the matrix element |(++|-ffint| )| in 

Eq. (|43)) . The eigenstates of the one-qubit part J2i=x y z ^^'^^ ^^'^ local Hamiltonian Hi^o are given by 



|+)=cos-|0)+e^^sin-|l), 
h)=-sin^|0)+e'*cos^|l). 



(48) 



where we define |0) and |1) as the eigenstates of and represent {h'^}i=x,y.z as {h sin 9 cos (p, h sin 9 sin (l),h cos 9} in 
the polar coordinate. We can thereby express the eigenstates |++) and | ) of -ffLO in the forms 



I++) = cos y COS y |00) + COS y siu y e"^^ |01) 

+ sin y COS y e^-^i 1 10) + sin y sin ^e*(^i+'^2) 1 11), 
I ) = sin y sin y |00) - sin y COS y e*^^ |01) 

_ COS ^ sin ^e^*MlO) + cos ^ cos ^e^(*i+*^) |11). 
22'' 22 ' ' 

We therefore have the matrix element (++|iJint| ) in the following form: 

{ + + \Hn.t\—) 

=Jz sin 6*1 sin 02 + (Jx - Jy) 



2^sin2^e-^(^^+*^) 
2 2 



-cos2^cos2^e*(^i+^^ 
2 2 



cos2^sin2^e^(*-*^) 
2 2 



.sin2^cos2^e^(-^^+^^) 
2 2 



(49) 



(50) 
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In the cases of {Jx, Jy} > > and > Jz > {Jx, Jy}, the upper bound of |(++|i?int| )| is given by 

|( + + |i/int|— >| <|J. + J,|. (51) 

We prove this inequahty in the cases oi > Jy > Jz > 0; we can prove the other cases in the same way. First, 
|(++|i?int| )| satisfies the following inequality: 

l( + +|i?inth- )| 

<\Jz \ sin 6*1 sin 6*2 
+ \Jx - ^ylfsin 

+ \Jx + Jy\^ 

\ T \ ■ a ■ a : \ T , 1 1 + cos 6li cos 6*2 
==1 sm6'i sin6'2 + \ Jx - Jy\ ^ 

1 — COS 9i cos 62 



2 ^ 
2 


sin^ 


f2 

2 " 


{- cos^ 


2 


cos^ 


02 
2 


2 ^1 

' y 


sin^ 


O2 
2 " 


f sin^ 


01 
2 


cos^ 


02 

2 



+ T ^ (52) 



By utilizing the fact that Jx ^ Jy ^ Jz ^ 0, the inequality (|52p reduces to 

|( + + |i?int| )| <>/2; - Jy COS COS ^2 + sin sin ^2 

<\Jx + Jy\- (53) 

The inequality (jSip becomes an equality when we choose {^1,^2,01,02} as {0, 7r,0, 0} for example, or in the Carte- 
sian coordinate {/if , /if , /if , /i^ , ''-i} = {0, 0, /i, 0, 0, — /i}. Then, the optimizing local parameters are given in the 
form of {/if , /i^, /if , /if, /12, ^2} = {0, 0, /lop, 0, 0, — /lop}. Moreover, Eqs. pT|) and can be given by substituting 

|(-|-+|iJint| )| with \ Jx + Jy\ in Eqs. ^ and (gS]). 

Finally, the leading order of Lambert's W function —W-i{—x) is log a; [2^. Therefore, the leading order of Eq. p6)) 
gives Eq. (HU, which then results in Eqs. ([T5|) and p^ . This completes the proof of Theorem 1. 



B. Optimization at arbitrary temperatures 

It is difficult to generalize Theorem 1 to arbitrary temperatures. However, we can present the following Theorem 2. 
Let us now parametrize the local fields as follows: 

{/if,/i?,/if,/if,/i^,/i^} 

={/if , /if, h^{i + e), /if , hi ~h^-{i - 0}, (54) 

or 

Theorem 2. When we express the negativity as a function of the local parameters 

{/if, h\, h^{l + ^), /if, /12, —h^{l — 0}, the following equation holds at arbitrary temperatures: 

dN _ dN _ dN _ dN _dN _ 
dJ^ ~ dh^ ~ dhf ~ dh^ " ~ 

at {/if, /if , /if , /i|, e, /i^} = {0, 0, 0, 0, 0, h}. (56) 

This theorem means that the form of the optimizing local parameters in the high-temperature limit, 
{/if , /if , /if , /if , /if , /if } = {0, 0, /lop, 0, 0, — /lop}, also gives an extremal value of the negativity at arbitrary temper- 
atures. 

Proof. To prove this theorem, we firstly calculate the perturbation of the negativity due to an infinitesimal variation 
of the local parameters at arbitrary temperatures. If it always vanishes, Eq. ([55)1 is proved. We first derive the 
perturbation of the density matrix due to an infinitesimal variation of the local parameters, from {0, 0, /i, 0, 0, —h} to 
{5hf, (5/if , /i(l -|- (5^), 5/if , (5/if , — /i(l — <50}- This means the perturbation of the form 

Htot = HZ+SHi^o, (57) 
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where 



i—x^y^z 

is the total Hamiltonian with the local parameters {0, 0, h, 0, 0, —h} and 

SHw = ('^^I'^i + '^^2^ ® ^2) 



(58) 



(59) 



is the infinitesimal variation of the local Hamiltonian. Equation (1211) gives the perturbation of the density matrix Sp 
as 



Sp 



-;3(H°p+5-Hlo) p-PH'Z 



Z + SZ 



z 



where pop — e l^^tat jZ and 



,5Z = tr(^-/3^' e-^(i-^'^°-<5HLoe-^"^°-dx^ . 



(60) 



(61) 



Then, the perturbation of the negativity in Eq. (|32|). (5A^ = — 2(0_|(5p-^i|0_), is given as 



<57V = _^iV(pop) 



2/3 



^^^^tr [J^ e-'5(i--)^-5HLoe-'^^^-dx 



tr 



tr 



^1 -/3(l-x)//° 



SHhoe 



dx 



LO 



dx 



where 



/3 



iV(pop)(/®/)+2(|0_)(0_|) 



(62) 



(63) 



and we used the identity N{pop) — —2{(j)-\p'^^\(j)-) as well as tr(A^ii?) — ti-{AB^^). We will prove that the integrand 
of Eq. dMI, 



tr 



always vanishes for {hf, h\, /if, /ij: ^z] = {0, 0, 0, 0, 0, K}. 
We prove in Appendix C that the operator 



(64) 



(65) 



has the same symmetry as the Hamiltonian H°^^ in Eq. (|58)) . and thereby must be expanded in terms of the Pauli 
matrices in the form 



(66) 
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where qoo, Qzo and qa are appropriate coefRcients. Therefore, we can calculate Eq. ((64|) to have the following equation: 

1 ■ 




qool «) / + gzo(o'i 0-2) 



i—x^y^z 



■i—x,y 



(67) 



A straightforward algebra, such as tr((Tf erf) = 0, yields that Eq. (|57)) vanishes. This means that the perturbation 
of the negativity due to the infinitesimal variation of the local parameters {Sh},., Shy, Sh'^, Shy, S£,} always vanishes at 
{/if, h\, /127 fT-z} = {O7 0, 0, 0, 0, h}. This completes the proof of Theorem 2. 

To extend Theorem 1 to arbitrary temperatures, we assume the following hypothesis: 
Hypothesis 1. The local parameters oj the form {/if, h\, /if, /i|, /if, /i|} ~ {0, 0, hop, 0, 0, —hop] give not only an extremal 
value but also the maximum value of entanglement at arbitrary temperatures. 

We numerically tested this hypothesis using determinant-based entanglement measure n{p) [2l[, which is given as 



7r(p) 




for > 0, 
for < 0. 



(68) 



Though this entanglement measure is not a full entanglement monotone, it provides tight lower and upper bounds for 
other entanglement measures including the negativity and the concurrence. In addition, det p^^ is expressed in the 
form of a polynomial and is much easier to maximize numerically than the concurrence and the negativity. Utilizing 
this measure, we tested Hypothesis 1 by numerical optimization for various kinds of interaction at various temperatures 
and found it always satisfied. In the following, we will assume Hypothesis 1 and conclude that {/if, /i^ , /if , /if , /i| , /ij } = 
{0, 0, /lop, 0, 0, —/lop} is the globally optimizing solution at any temperatures. 

For the local parameters {/if, /if, /if , /if , /i|, /12} = {0, 0, h, 0, 0, —h}, the density matrix Zp^^ is given at arbitrary 
temperatures in the basis of the eigenstates of erf (g) cr|, {|00), |01), |10), |11)}, as 



ai 








bi 


-62 


63 





63 


61 + 62 


a2 









(69) 



where 



ai 

61 
63 



cosh;3Ji, 02 = - 



s''^^(J^ + Jy)sinh/3J2 



cosh 13 J2, 62 



J2 

2/ie'^'^= sinh ^ J2 
J2 ' 



-er^'^- sinh/3Ji, 
Ji = |J,-Jj, J2 = y'4/i2 + (J, + Jy)2 



(70) 



Its eigenvalues are 



{ai - \a2\,ai + |a2|,6i + ^bl + bl,bi - ^bl + bl^. (71) 
In Appendix D, we will prove that only ai — |a2| can have a negative value for { J^^, J^} > Jr > and > Jz > {Jx, Jy}- 



Therefore, the optimized negativity is given by 

N{J^,Jy,Jz,h,l3) = max(7V,0), 



(72) 
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where 



N 



Z 



2 »i ~ |q2| 
z 

e-^-'-' cosh/3Ji - {el^-'-lJx + Jy\ sinli /JJaj/Js 
e^'^-'' cosh/3Ji + e'^-'' cosh/3J2 
2e~'^'^^ cosh^Ji + 26^-^^ cosh/3J2. 



(73) 



We find from tliis expression that we can always make the negativity positive by clroosing an appropriate vakie of h. 

The remaining task is to find the value of the optimizing field hop at each temperature. We will do it analytically 
in the low-temperature limit /3 — )■ oo in Sec. III.C as well as do it numerically rigorously for a wide range of the 
temperature in Sec. IV. 



We now discuss the optimization problem in the low-temperature limit. 

Theorem 3. In the low-temperature limit /? — >■ oo, the optimized entanglement approaches to 1. The optimizing 
parameter hop approaches to when we choose the optimizing parameters as {0, 0, /lop, 0, 0, — /lop}. 

Proof. We need to consider the three cases, namely the cases where the ground state of i?int is non-degenerate, 
doubly degenerate and triply degenerate. The eigenvalues and the corresponding eigenstatcs {|V'i)}f=i of ^int 

are given by the following: 



As has been described in Sec. II, we consider only the cases of { J^-, Jy} > Jz > and > J2 > {Jx, Jy}- 

In each case of {Jx, Jy} > Jz > or > J^ > {Jx, Jy}, the ground state of Hint is non-degenerate, and ei or £2 
is the ground-state eigenvalue, respectively. In these cases, the ground state is a Bell state and it is clear that its 
entanglement is maximum. In other words, there is no need to optimize it further and — 0. We will see in Sec. IV 
that, in this non-degenerate case, there is indeed a finite range of the temperature where the negativity is maximized 



In each case oi > Jz = Jx > Jy and > Jz ~ Jy > Jx, the ground state of Hmt is doubly degenerate and £2 = £4 
or £2 = £3 is the ground-state eigenvalue, respectively. In the case > Jz = Jx = Jy, the ground state of Hmt is triply 
degenerate and £2 = £3 = £4 is the ground-state eigenvalue. In these degenerate cases, the ground states are mixed 
states and their entanglement always vanish. However, we can resolve the degeneracy of the ground states by an 
infinitesimal local Hamiltonian. 

We hence employ Hypothesis 1 and put {h^,h\,h\,h2,h\,h2} — {0, 0, ft-op, 0, 0, — /lop}- We then calculate the 
asymptotic behavior of the optimized entanglement in the low-temperature limit /3 — 00. Below we will derive 



C. Optimization in the low-temperature limit 




(74) 



for /?°P = 0. 




1 + log 2/3 J 
PJ 



as /3 — 00. 



(75) 
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in the doubly degenerate cases, where we defined J = | Jj; + Jj,|, and 



~ y ^ log 4/3 J as /3 ^ oo 

^„p^l_i±i^ as P^oo. (76) 

in the triply degenerate case. In both cases the optimizing parameter hop is infinitesimal and the optimized negativity 
iVop approaches to 1 in the low-temperature limit ^ — >■ oo, although the forms of hop and A^op are slightly different in 
the two cases. We will see in Sec. IV that, in these degenerate cases, there is indeed no finite range of the temperature 
where the negativity is maximized without local fields. In other words, we need a non-zero value of hop at any non-zero 
temperatures. 

Now we derive Eqs. ([75]) and ([75)1 . We start from Eq. ([75]) under Hypothesis 1. In the doubly degenerate cases 
> Jz = Jx > Jy and > Jr = Jj, > Ja;, we can approximate Eq. ([75]) as 

^ ^ _e^-^^) - (e'3(^^+-'^)| J. + Jy\)/J2 



g-/3(,/,-,/i) ^ g/J(J, + J2) 



I _j_ g/3(2J, + J2-Ji) 



(77) 



in the low-temperature limit /3 — >■ oo, where we used the facts that 2 cosh/3 Ji ~ e'^-^^ , 2 sinh/3J2 ~ e'^'^^ and 
2cosh/3J2 ~ e^''^. Moreover, in these doubly degenerate cases, 2Jz + J2 — Ji is either 2Jz + J2 — Jx + Jy or 
2J^ + J2 + Jx ^ Jy, which are summarized to J2 — | J^; + Jy|. Then Eq. (|77p reduces to 

~ ^ -l + e^^J/{X + J) 

^- TT^ ' (^^^ 

where 



X = J2 - J = •y/4/i2 + (J, + j^)2 _ I + (79) 

We first prove that X — )- and /3X — > 00 is a necessary and sufficient condition for /V — > 1 in the low-temperature 
limit P 00. In order to prove this, we calculate the value of 1 — iV as follows: 

1 + ef^ 
2e'0^ + X/{X + J) 

= ' 

Because X > and < e"''^ < 1, we have X/{X + J) > and 1 < 1 + e'^^^ < 2. Therefore, the necessary and 
sufficient condition for 1 — /V —> in the low-temperature limit is 

/3X 00 and X ~> as p ^ 00. (81) 

In such cases, the negativity can be maximized to 1 in the low temperature limit /3 — > 00. 

Let us now calculate the optimizing parameter Xop. From the extremal condition for Eq. (|78p . 



dN e''^(/3X2+3/3JX + 2/3j2- J~e'3^") ^ , . 

— U, (»/) 



dX {l + ef^)^J + Xy 

we obtain 



pXop = log(— + ipXop + 2/3 J - 1) . (83) 
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Because of the condition (|81l) . Eq. ((83|) reduces to 



-log 2/3 J (84) 

in the limit /3 — >■ oo. We thus have 

X.-^^. (85) 

which indeed satisfies ()8ip . The optimizing parameter /lop is thereby obtained in the form 

1 



/lop = 2 V -'^op + 2 J ^op 



^y^log2/3J, (86) 
where we utilized Eq. (|79p to derive the first equality. Moreover, the optimized negativity is given by 



iVop 



e-/3Jfop + 1 
^ (l-e-^^°p)(^-e-^^- + l-^ 

~ 1 - ^ - 2e-^^- 
^^_l + logW 

where we used Eq. ([5T|) upon moving from the first line to the second line. Thus Eq. ([75]) is proved. 

In the triply degenerate case > Jz = Jx = Jy, we have Ji = 0, and thereby we can approximate Eq. (j73p as 

^^_2e-^./._(e^(./.+./.)|J, + J^|)/J, 



2+(e/3(2./.+J.)|j^ + jj)/j^ 

2 + e/3(2,/, + ,72) 



(88) 



in the low-temperature limit /3 — > oo. where we used the facts that cosh/3Ji ~ 1, 2 sinh /3 J2 ~ e^"^^ and 
2cosh/3J2 — e^''^. Moreover, in this case, 2Jz + J2 is equal to J2 — \ Jx + Jy\, and therefore Eq. (|55| reduces to 

^ = -2 + e^"j7(^ + J) 
2 + ' 

where X and J arc defined in Eq. ([75)1 . From the extremal condition dN/dX — 0, we obtain 

^2 



Xop/? - log(— + 6/3Xop + 4/3 J - 2) 



J 

log 4/3 J, (90) 



where we used the same logic as the one with which we derived Eq. ((84)) in the doubly degenerate case. In this way, 
the optimizing parameter hop and the optimized negativity N^p are given as 



hop^J^\ogAf3J, (91) 
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and 



.l-i±i^^. (92) 
I3J 



Thus Eq. ([75]) is proved. This completes the proof of Theorem 3. 



D. Negativity and Concurrence 

We here mention the relationship between the negativity and the concurrence [l^. The concurrence is also an 
important entanglement measure. Concerning the optimization problem of the concurrence, we can only prove that 
the negativity TV and the concurrence C have the same value for the local parameters {0,0, h, 0,0, —h} with an 
arbitrary value of h; namely, 

N{0, 0, h, 0, 0, -h, P) = C(0, 0, h, 0, 0, -h, 13). (93) 

This equation is proven by the theorem in Ref. [30j , which says that the concurrence is equal to the negativity iff the 
eigenvector of p^^ corresponding to its negative eigenvalue is a Bell state up to local unitary transformations. 

For the local parameters {0,0,h,0,0,~K], the density matrix Zp^^ is given in Eq. and only the eigenvalue 
N = ai — \a2\ can be negative. For {J^, Jy} > Jz > and > Jz > {Jx, Jy} the eigenvectors of Zp^^ corresponding 
to the eigenvalue TV = oi — |a2| is (|00) + |11))/V2 and (|00) - |11))/V2, respectively, both being a Bell state. In 
the case of TV > 0, the concurrence must be equal to the negativity because the eigenvector of p'^^ corresponding to 
its negative eigenvalue is a Bell state. In the case of TV < 0, the negativity TV = max(TV, 0) is equal to and the 
entanglement does not exist. Therefore, the concurrence and the negativity are both equal to 0. This completes the 
proof of Eq. ^ 



IV. HIGH- AND LOW-TEMPERATURE PHASES 



In the present section, we calculate the optimizing local Hamiltonian and the optimized entanglement numerically 
rigorously. After the analysis in Sec. Ill, we here set {h\, h\, /if , /if , /i|, /i|} = {0, 0, /lop, 0, 0, —hop}- In the calculations 
below, we will see that there are two kinds of temperature range, which we refer to as the high- and low-temperature 
phases. We will find that in the low-temperature phase the optimizing local parameter /lop vanishes, whereas in the 
high-temperature phase it dose not. We start from Eq. ([75]) with the optimizing parameters {/if, h\, h\, /if, /if, /i|} = 
{0, 0, /lop, 0,0, —/lop}. The parameter /lop which maximizes the negativity can be calculated from 



1 dN 
h'dh 



oi\Jx + Jy\{PJ2 cosh/3J2 — sinh/3J2) 

h=h„„ 



+ J2 13 smh 13 J 2 7— -J2/3 H 

cosh/3Ji v 2 

=0, (94) 

where the factor l//i is added to remove the trivial solution of h = 0. In Fig.[Tl we show the optimizing local parameter 
/lop in the cases of { J^, Jy, JJ = {1/3, 1/3, 1/3}, {1/2, 1/3, 1/6} and {-1/2, -1/4, -1/4}. See Appendix A for the 
convergence of /lop to the asymptotes ([TT]) and p3|) . 

In the high-temperature phase, Eq. (j94p has a non-trivial solution of hop > 0, while in the low-temperature phase, 
Eq. (|94p has no solutions and the optimizing value hop is zero, which is the trivial solution of dN/dh = 0. Therefore, 
the boundary temperature Tc between the high- and low-temperature phases is a solution of 

l dN{Jx,Jy,J.,h,P) 

hm ^ = 0. 95) 

h^o h dh 

The boundary temperature is defined for each interaction Hamiltonian -f/int- 

In Fig. [21 we show the boundary temperature in the cases of { J^,, Jy} > > and > Jz > {Jxt Jy}j which 
correspond to all kinds of interaction thanks to Lemma 1. We calculated the data in Fig. [5] from normalizing 
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FIG. 1: Numerically rigorous solution of the optimizing local parameter hop', the solid line for {Jx, Jy, Jz} = {1/3, 1/3, 1/3}, 
the dashed line for {J^, Jy, J^} = {1/2, 1/3, 1/6} and the thick line for {J^, Jy, J^} = {-1/2, -1/4, -1/4}. The number of data 
points is 3000 for each case. The boundary temperatures Tc between the high- and low-temperature phases are 0.8168 • ■ ■ , 
0.6803 • • • and for {1/3, 1/3, 1/3}, {1/2, 1/3, 1/6} and {-1/2, -1/4, -1/4}, respectively. 




FIG. 2: (color online) The boundary temperature Tc between the high- and low-temperature phases, (a) for the antiferro- 
magnetic case {Jx, Jy} > Jz > and (b) for the ferromagnetic case > J, > {Jx, Jy}- The point of origin is (1, 1), which 
corresponds to the isotropic Heisenberg interaction. In (a), the boundary temperatures are 0.8168 • • • , 0.1292 ■ ■ ■ , 0.1292 • ■ ■ , 
0.5735- ■■ and 0.6208 •• • at (1,1), (1,10), (10,1), (10,10) and (5,5), respectively. The maximum temperature is 0.8168 •• • 
at (1, 1), which is the XXX point. In (b), the boundary temperatures are 0, 0, 0, 0.4126 ■ • ■ and 0.3188 ■ • ■ at (1, 1), (1, 10), 
(10, 1), (10, 10) and (5, 5), respectively. The maximum temperature is 0.5184 • • • at lima;_>oo {x,x), which is the XX point. 

the interaction parameters so that ||i?int||2 — 1, where || II2 is the spectral norm. From Fig. [21 we see the following 
properties. First, the boundary temperatures Tc are higher in the antiferromagnetic cases {Jx, Jy} > Jz > than in 
the ferromagnetic cases > > {Jx-, Jy}- Second, in the antiferromagnetic systems, the boundary temperature Tc is 
maximum of 0.8168 • • • for the isotropic Heisenberg interaction (the XXX model). Next, the boundary temperature 
Tc is zero in the cases oi > Jz = Jx > Jy and > Jz = Jy > Jx a.s well as the case of the ferromagnetic isotropic 
Heisenberg model, which means that the low-temperature phase shrinks to the zero temperature in these doubly 
and triply degenerate cases analyzed in Sec. HI.C. We have revealed in Sec HI.C that in the low-temperature limit 
/3 — > 00 the negativity is strictly 1 with no local Hamiltonian in the non-degenerate cases. The present calculation 
indeed shows that the low-temperature phase extends to a finite temperature in the non-degenerate cases. In the 
antiferromagnetic system, on the other hand, the boundary temperature is zero only in the case of the Ising model, 

Jx = Jz = 01 Jy = Jz = 0. 

Next, we consider the singularity at the boundary between the high- and low-temperature phases. In Fig. |3l we 
show the optimized negativity, its first derivative and the purity tr(yO^) in the case of {Jx, Jy, Jz} = (1/2, 1/3, 1/6}. 
We also consider the entanglement enhancement, which is defined as the difference of the entanglement between 
the optimized entanglement and the entanglement under no local Hamiltonian, namely N{H'^q) — N{Hi^o — 0)- We 
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FIG. 3: Plots of (a) the negativity; (b) its first derivatives; (c) tlie purity, for J^c — 1/2, Jy — 1/3, Jz = 1/6. The solid line is the 
optimized entanglement and the dashed line is the entanglement enhancement defined in the text. The boundary temperature 
is Tc — 0.6803 ■ ■ ■. At T = 1.185 • • •, the entanglement enhancement is maximum, where the value is 0.1480 • • • . In (b), we 
obtained the data points by the finite-difference method. In (c), the solid line is the purity of the optimized state and the 
dashed line is the one under no local Hamiltonian. The minimum value of the purity is 0.5087 • ■ ■ at Tc. 



numerically rigorously calculated the data in Fig.[31[a) using ([M)) . and the derivatives by the finite-difference method. 
Figure [D^b) shows that the second derivative of the negativity is not continuous at the boundary and Fig. ^c) shows 
that the first derivative of the purity is not continuous at the boundary. On the other hand, there is no singularity at 
the point of T = 1.185 • • •, where the derivative of the entanglement enhancement is not continuous. 

The emergence of the high- and low-tempcrature phases is due to the following reason. First, the entanglement 
enhancement by addition of the local Hamiltonian comes from the fact that a local Hamiltonian increases the purity 
and suppresses the entanglement loss caused by thermal mixing, which is demonstrated in Fig. ^c) . On the other 
hand, too strong magnetic fields bring the quantum system close to a direct product states and hence destroy the 
entanglement. These two effects compete to give rise to the two phases. In the low-temperature phase, we do not 
need a magnetic field because the purity is already high. In the high-temperature phase, on the other hand, we need 
a magnetic field because the thermal fiuctuation decreases the purity. The transition from the low-temperature phase 
to the high-tcmperature phase means that the enhancement of the entanglement due to the increase of the purity 
becomes predominant compared with the entanglement decay caused by the magnetic decoupling. 



V. SUMMARY AND CONCLUSION 



We have analytically and numerically rigorously studied thermal states of quantum systems where two qubits in- 
teract under a local Hamiltonian Hi^o f^nd have determined the local Hamiltonian Hi^o which maximizes the thermal 
entanglement under a fixed interaction. As a result, we have found that the interaction Hamiltonian can be transformed 
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into the Xy ^-exchange interactions whose parameters are either antiferromagnetic as { J^;, Jj,} > > or ferromag- 
netic as > Jz > {Jx, Jy} and that the optimizing local Hamiltonian always takes the form of hop{af (g) / — J (SD (t|), 
where hop depends on the temperature. In addition, we have proved that the optimized entanglement docs not vanish 
at any temperatures and decays slowly according to l/(nogr) at high temperatures. We have also found that in 
the low-temperature phase the entanglement is maximum without any local Hamiltonian and have investigated the 
interaction dependence of the boundary temperature of this range. Indeed, the low-temperature phase shrinks to 
the zero temperature point if the interaction Hamiltonian has degeneracy. At the same time, we have discovered a 
singularity of the optimized entanglement at the boundary temperature, where the second derivative is discontinuous. 

In conclusion, our work has revealed general properties of the thermal entanglement of interacting two qubits, 
though we have assumed a numerically confirmed hypothesis. The concept of high- and low-temperature phases is 
an interesting property in that it is based on the response to external manipulation of local Hamiltonians. It is likely 
that we can find more interesting properties of entanglement in this regard. In future, we plan to investigate two 
qubits which interact indirectly or general bipartite systems. 
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Appendix A: Numerical comparison 

In the present Appendix, we compare the asymptotes in Eqs. (fTTj) and ((T3|) with those in Eqs. (fTSj) and (|16p in the 
case of {Jx, Jy, Jz} = {1/3, 1/3, 1/3}. In this case, Eqs. (fTTj). ([13]), (fTS]) and (|T6| . respectively, reduce to 



/3 



(Al) 



"El - (A3) 

m = ^si^. (AD 

where the subscripts denote the equation number of the corresponding asymptotes. In Fig. 21 we show the comparison 
of these asymptotes with the numerically rigorous estimates of hop and A'op obtained from Eq. (|94| . We can see that 
the convergences of ^ fY5l /^°p ^^'^ A^ jYSl /A^pp are very slow, while the convergences of /iJXTj/^op ^^'^ iV jY^ /jVQp are 
much faster. The convergence of N{h^YT^)/Nop, where N{h) is given in Eqs. ([7^ and ([75)1 . is even faster than that of 
N^/Nop; at T = 100, the values of iV(/i|ii})/iVop and N^/Nop are 0.999998 and 0.9994, respectively. 



Appendix B: Lemma 2 in degenerate cases 



In the proof of Lemma 2, we left out the cases of C = or = ±1 in Eq. ((23)) . In the present Appendix, we prove 
that Lemma 2 still holds in these cases. First, the general form of Eq. (PT|) is given in the basis of as follows: 



Zp 



012/12^ 031/31/3 a32/32/3\ 



021/21/3 e 



2Ch' 



041/41 /3 042/42 /3 



013/13/3 oi4/i4;5 e 2'^''' 034/34/3 
\023/23/3 044/42/3 a^^zfiiP e-^''' / 



(Bl) 



where {oy} are constants of order 1 and {fi_iu} are defined in Eq. (P^ . Note that on the diagonal of Eq. (|B1[) . the 
second term of Eq. (PT|) is neglected in comparison to the first term. In the cases of C 7^ and C 7^ ±1, Ec^. (jBl[) 
reduces to Eq. (^5)) . 
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h{ii)/hop > h{i5)lho 
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FIG. 4: The comparison between the asymptotes in Eqs. (Illf) and (|13|l and those in Eqs. (|15|l and (|16|l . (a) for the ratios 
/i jYY] //iop (solid line) and ftjXgj/^op (dashed line), where /iJYT) ^'^'^ '^(TSj derived from Eqs. (|AH) and (|A3[) . respectively, 

and /lop is the numerically rigorous value calculated from Eq. (|94|) . At T = 100, the values of /li^ /hop and /iJ^^/Ziop are 
1.0007 and 0.4438, respectively, (b) for the ratios N ^^ /Npp (solid line) and -^ |T6] /-^qp (dashed line), where A^ |]"3] and Afjjgj 
are derived from Eqs. HA2p and (|A4|I . respectively, and A^op is the numerically rigorous value calculated from Eq. (|94|l . At 
T = 100, the values of N^/Nop and Ajjgj/A'op are 0.9994 and 2.494, respectively. 



In the case of C = 0, we have {E[, E!^, E^, E!J = {2/i', 0, 0, -2/i'} and /23 = fsi = e"^^ = 1, and hence Eq. ((BT|) 
reduces to 



/ ^2h' 






032/3 




1 


041 ^, 


042;^ 






1 


034-^ 


\ 023^ 


024 77 


P 

«43 7r7 


g~2h' 



(B2) 



In this case, the product of the diagonal elements (PD) of Zp^^ is 1, whereas the maximum of the absolute val- 
ues of the products including ofi-diagonal elements (POD) is of order e^^ /J^/Zi'^, which comes from the product 
,§e — ^ ^-2h ^ Therefore, it is necessary for detp^^ < that the order of e^^ fi"^ /h'"^ is greater 



^2h' ^ „ Pe" 

-e X 041 — 



X ai4^ 



or of order 1, which leads to 



logl//3 



(B3) 



as in Eq. (P^ . Thus, Lemma 2 is proved in the case of C = 0. 

The proofs for the cases of C = 1 and ^ = —1, or the cases of {/ii,/i2} = {2/i,0} and {/ii,/i2} = {0,2/i}, are es- 
sentially the same. We here present the proof only for the case of C = 1- In this case, we have {E[, E2, E'^, E'^} = 
{2h', 2h', -2h', -2/1'}, /12 = /21 = e2''' and /34 = /43 = e^^'*', and hence Eq. dE]) reduces to 



2h' „2h' 



2h' 



a2i/3e 

/3e- 



2h' 



ai3- V 

\a23V 



e 
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041 V- 

^-2h' 



032^ \ 
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043Pe 
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(B4) 



The PD of Zp^^ is 1, whereas the maximum of the absolute values of the PODs is of order e^^ or of order 



-e^^' X 041 ^^^7— X 014 ^^p— X e and 032 



Pe' 



X 041^^ X ai4 



Pe' 



X 023^^, re- 



g8/i pi jy^ii^ which come from 

spectively. Therefore, it is also necessary for det p^^ < that e'^'' 0^ /h''^ is greater or of order 1, which again leads to 
Eq. (|B3|) . Thus, Lemma 2 is also proved in the case of C = 1. 
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Appendix C: Proof of Eq. (f66|) 

In order to prove Eq. (|66|) . we begin with the standard operator expansion of an arbitrary 2® 2 operator Q: 

Q = \ E (CI) 

i.j—Q.x.y.z 

where cr^ = CT2 = / is the two-dimensional identity operator. The coefficients qij are given by 

=tr(ga* (C2) 

because tr(/ ® I) — A and the other terms are traceless. 

Symmetries that the Hamiltonian (j58p possesses ehminate many of the coefficients {qij} of the expansion of operators 
with the same symmetries, such as exp(— /3a;i/°^t). First, a straightforward calculation shows that the Hamiltonian ([55)1 
commutes with the global phase flip 

{/flip = e''(''/2)'"i ® e'('^/2)'^2 

= -al®(jl. (C3) 

This operator flips the signs of and u^. For an operator Q that commutes with C/gip, the coefficients 
{qox, qov,qxO, qyo, qxz,qyz,qzx,qzy} vanish. For example, we have 

qxz = tr[Q{a=^ ® al)] 

= tr[C/flip"^QC/flipC/flip"H(T? ® o-2)%ip] 
= tr[Q((-aJ=) ® 

= -qxz = 0. (C4) 

The same argument gives qox = qoy = qxo = qyo = qxz = qyz = qzx = qzy = 0. 

Next, the Hamiltonian (j58p is a real matrix in the basis. Noting that only has imaginary elements in this 
representation, we have, for an operator Q with the symmetry Q* = Q, 

(9.,)*=tr[Q*(Kr®(^!)*)] 

= ~qxy (C5) 
On the other hand, the Hermiticity of an operator Q is followed by 

(g.,r =tr[(K)t®(a|)t)Qt] 
= tr[Q{af ® ay)] 

= qxy (C6) 

The above argument shows q^y = qyx = 0. 

Finally, the Hamiltonian ([55)) is symmetric with respect to the following set of operations: 

= -K®^f)Pi2, (C7) 

where P12 is the permutation of the spins 1 and 2. The operator af (g) erf flips the signs of erf and erf but the 
permutation P12 makes the signs back to the original ones, because the local fields are in the opposite directions in 
the Hamiltonian ([58)) . For an operator Q that commutes with C/12, we have 

g.O =tr[Q(crf «)/)] 

= tr[[/i2-ig;7i2C/i2"'K ® i)Ui2] 

= tv[Q{I®{-cjm 

= -go.. (C8) 
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To summarize, an operator with the same symmetries as the Hamiltonian (|58p is expanded in the form 



E 



(C9) 



In the operators e ^'-^^tot and e '^^^ ^'^to't have the same symmetries as the Hamihonian H°^^ and hence are 

given in the form (jC9l) . 

Since the density operator p = e~^^*°* is given in the form (jC9p . the partial transpose p^^ is also of the form (|C9[) : in 
the (T^ basis, the partial transpose Ti only flips the sign of af and hence changes only the sign of Qyy in the expansion, 
not the symmetries nor the form of the expansion. 

The state |(^_) is a non-degenerate eigenstate of the operator p-^^ if the minimum eigenvalue A_ is negative. Suppose 
that the operator p'^^ commutes with a symmetry operator U. Then the projection operator should have 

the same symmetry. This is shown as follows. Since we have 



(CIO) 



and |(/)_) is non-degenerate, the vector U\(j)-) must be the same vector as \(f)^) except for a phase: U\(j)-) = e'^|0_ 
Therefore, the projection operator commutes with U if the negativity is non-zero. This means that 

as well as have the same symmetries as the Hamiltonian H°^^ and are expanded in the form (|C9p . 

We thereby arrive at the conclusion that the operator 



iV(Pop)(/®/)+2(|,^_)(0_|) 



Ti 



(l-^)fftTt 



(Cll) 



has the same symmetries as the Hamiltonian H^^^ and hence is expanded in the form (jC9 



Appendix D: The Eigenvalues of (|69|) 

In this section, we prove that in the eigenvalues of the matrix (|69p . only oi — \a2\ can have a negative value 
for {Jx, Jy} > Jz> and > Jz > {Jx, Jy}- The four eigenvalues are given in ([7T|) . Because ai > 0, \a2\ > and 
b2 + > 0, obviously have 



ai + |a2|>0, bi + Jbl + bl>0. (Dl) 



Therefore, we only have to prove that hi — + > 0- 

First, we prove this inequality for h ~ 0. For h ~ 0, the eigenvalue bi — 1/62 + ^3 reduces to 



bi-^bi + q 

= e^-^' cosh[/3(J, + Jy)] - e-^-^' sinh[/3| J, - Jy\] 

_ g^P{\J.-Jy\-J.) f^Pi-\J.-Jy\-J.)y (132) 

For {Jx, Jy} > > 0, we have 

gP{J.+Jy+J.) „ g/3(|Jx-J„|-J.) > 0^ (D3) 



which leads to 61 — 1/62 + 63 > 0. For > Jz > {Jx, Jy}, we have 



gPi-J.-Jy+J.) _ ^P{\J.-Jy\-J.) 



'2e-/3'^« smh[l3{-J^ + Jz)] for > Jz > Jx > Jy, 
2e-P-^^ sinh[/3(- Jj, + Jz)] for 0> Jz> Jy> Jx- 



(D4) 
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Because —Jx + Jz > and — J.„ + > 0. 



(D5) 



ioT > Jz > {Jx, Jy}, which also leads to &i 



Next, we prove bi 



hi 



h\ > 0. Thus, hx - \Jh\ + h\ > is proved for h 



0. 



The value of h\ — h\ — h\ is calculated as follows: 

}? — }? — }? 
Oi "2 "3 



6| > for arbitrary /i, which is equivalent to 61 — \Jh\ + h\ > because 61 + \/&2 + H > 0. 



= 2/3 



■^^ ^cosh^ /3J2 - ^ sinh^ PJ2^ - e-^^-^' sinh^ /3 Ji 



1+^1 



sinh^ /3J2 
^sinh2/3J2 



'2,3,/^ 



sinh /3Ji 



~2/3J^ 



•^2 



sinh /?Ji 



(D6) 



Only the second term depends on h through J2 = V4^?TTXTJyF- The term (sinh/3J2/J2)2 is a monotonically 
increasing function of J2 for J2 > 0, while J2 is a monotonically increasing function of . Therefore, &f 
also a monotonically increasing function of h^. Since we already proved that 



obtain 6f 



63 > for any values of and thus 61 — + 6| > is proved. 



63 is positive for h 



0, we 
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